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Abstract. We study the validity of inspiral templates in gravitational wave data 
analysis with Advanced LIGO sensitivity for low mass binary black holes with total 
masses of M < SOMq. We mainly focus on the nonspinning system. As our complete 
inspiral-merger-ringdown waveform model we assume the phenomenological 

model, “PhenomA”, and define our inspiral template model (Xmerg) by taking the 
inspiral part into account from IM1Z up to the merger frequency (/merg)- We first 
calculate the true statistical uncertainties using XMIZ signals and XMIZ templates. 
Next, using XMIZ signals and Xmerg templates, we calculate fitting factors and sys¬ 
tematic biases, and compare the biases with the true statistical uncertainties. We find 
that the valid criteria of the bank of X^erg templates are obtained as Merit ^ 24M0 
for detection (if M > Merit, the fitting factor is smaller than 0.97), and Merit ^ 26Mq 
for parameter estimation (if M > Merit, the systematic bias is larger than the true 
statistical uncertainty where the signal to noise ratio is 20), respectively. In order to 
see the dependence on the cutoff frequency of the inspiral waveforms, we define an¬ 
other inspiral model Xisco which is terminated at the innermost-stable-circular-orbit 
frequency (/isco < /merg)- We find that the valid criteria of the bank of X[sco templates 
are obtained as Merit ^ ISMq and ^ ITMq for detection and parameter estimation, 
respectively. We investigate the statistical uncertainties for the inspiral template mod¬ 
els considering various signal to noise ratios, and compare those to the true statistical 
uncertainties. We also consider the aligned-spinning system with fixed mass ratio 
(mi/m 2 = 3) and spin (x = 0.5) by employing the recent phenomenological model, 
“PhenomC”. In this case, we find that the true statistical uncertainties can be much 
larger than those for the nonspinning system due to the mass-spin degeneracy. For 
inspiral PhenomC templates truncated at /merg, the fitting factors can be better but 
the biases are found to be much larger compared to those for the nonspinning system. 
In particular, we find significantly asymmetric shapes of the three-dimensional overlaps 
including bimodal distributions. 
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1. Introduction 

The next generation gravitational wave detectors, such as Advanced LIGO [I] and 
Virgo, [2] are likely to allow us to observe the real signals in coming years |3j . Coalescing 
binary black holes (BBHs) are among the most promising sources of gravitational wave 
transients for the ground-based detectors jl]. A coalescing BBH system suffers three 
phases: inspiral-merger-ringdown (IMR). In the inspiral phase, the two compact objects 
move in quasicircular orbit mutually approaching driven by radiation reaction. In the 
merger-ringdown (MR) phases, the system reaches the ultra-relativistic regime, the two 
bodies merge to form a single excited Kerr BH and eventually that settles down into a 
Kerr BH. The gravitational waveforms from the early inspiral phase can be accurately 
obtained by the post-Newtonian (PN) approximation (refer to Ref. [S] for an overview on 
various PN approximants), and in the ultra-relativistic regime the accurate waveforms 
can be calculated by the numerical relativity (NR) simulations [Hj. Of course, from 
the NR simulations one can extract the complete IMR waveforms (e.g. see Ref. 0). 
However, since performing a long NR simulation is computationally very expensive, 
most of those simulations have been done only in the last few orbits. On the other hand, 
efforts to establish the analytic IMR waveform models for nonspinning BBHs have been 
made by several authors |Hl El [El E] by means of the hybrid IMR waveforms, which 
can be obtained by combining the PN inspiral waveforms [13] and the numerical 
MR waveforms milEllElinilEllE!. Further studies now allow us to have the Fourier 
domain spinning IMR waveform models for aligned-spin [201 El] and processing [22] 
BBHs. 

For low mass compact binaries, whose components consist of a stellar mass BH 
and/or neutron star, the inspiral phase is likely to have most of the signal power and 
accurate inspiral waveforms can be computed by the PN approximants. Therefore, 
the PN inspiral waveforms have been generally used in the ground-based gravitational 
wave search [23] and parameter estimation [211 m\. Especially, the stationary-phase 
approximated PN waveform model (called “TaylorF2” [S]) has been mainly used because 
the waveform can be given by an analytic function in the Fourier domain. However, if 
the MR phases are nonnegligible in the total signal power, the complete IMR waveforms 
should be used as templates in searches not to lose a fraction of the signal to noise ratio 
(SNR). As the binary mass increases, the contribution level of MR phases to the SNR 
tends to increase. For detection purposes, thus, one can choose a critical mass (Merit) 
considering both the computational advantage and the detection efhciency For the 
signals with masses M < Merit simple inspiral-only templates can be used to lower the 
computational cost although there can be a small loss of the SNR. While, for the signals 
with M > Merit complete IMR templates should be used not to lower the detection 
efficiency below a certain value (typically 90%). Buonanno et al [5] and Brown et 
al [2S] showed that Merit ~ 12 Mq for various PN inspiral template models with IMR 
signals computed by the EOBNR model. Ajith da also showed that Merit ^ 15Mq for 
3.5PN TaylorTl inspiral templates with phenomenological IMR signals. The latter is 
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consistent with one of our results, which we will see in section 

In parameter estimation, inspiral waveform templates can produce systematic biases 
due to the difference between two models for the IMR signals and the inspiral templates. 
Farr et al 123 investigated the overall trend of fitting factors and biases using EOBNR 
signals with total masses up to ~ 300 Mq and the TaylorF2 and the FOB inspiral 
template models by performing Monte Carlo simulations. A similar work was carried 
out by Bose et al |2H] using phenomenological IMR signals with masses in mi ^2 G 
[13,104] Mq and the 3.5 PN TaylorTl inspiral templates. Bose et al |28] calculated the 
overlap surfaces numerically to obtain the overall trend of fitting factor and biases, and 
compared the biases to those obtained by using the analytic approximation described 
in Ref. j29] in the low mass region mi ^2 G [5,20]Mq. In this work, we also investigate 
systematic biases and fitting factors of inspiral template models for the nonspinning 
system. We employ the phenomenological waveform model, “PhenomA” j9] as our 
complete IMR model, and consider inspiral templates constructed by taking only the 
inspiral parts into account from the original IMR waveforms. However, we focus on low 
mass BBHs in the range of M < 30Mq and access the validity of inspiral templates 
for both detection and parameter estimation in detail, giving the quantitative results. 
Thus, our work can complement the previous studies and present some interesting new 
findings, which we will discuss in section [3^ 

On the other hand, astrophysical BHs are likely to have spins, that can significantly 
affect both the search and the parameter estimation in GW data analysis. Therefore, 
our results for the nonspinning binary system cannot be directly applicable to spinning 
systems. In order to see the impact of the spin on our results, we consider the aligned- 
spinning system with fixed mass ratio (mi/m 2 = 3) and spin (y = 0.5) by employing the 
recent phenomenological model, “PhenomC”, and show a brief comparison between the 
nonspinning and the aligned-spinning systems. We find that the degeneracy between 
the mass ratio and the spin can significantly affect the overlap distributions; thus, the 
statistical uncertainties and the biases can be much larger than those for nonspinning 
binaries. 

This paper is organized as follows. In section]^ we briefly review PhenomA model, 
and describe how to calculate statistical uncertainties from overlap surfaces in the 
Bayesian parameter estimation framework. In section using the IMR waveforms, 
we investigate statistical uncertainties for signals with a SNR of 20. Next, assuming 
a bank of inspiral templates, we calculate fitting factors and biases, and investigate 
the valid criteria of the inspiral template bank for the efficiencies in detection and 
parameter estimation. A systematic study on the impact of the SNR on our analysis is 
also discussed. Finally, we present some results for aligned-spinning binaries to see the 
spin effect. A summary of this work is given in section with future works. 
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2. Phenomenological Waveforms and parameter estimation 


2.1. Phenomenological Waveforms 


In the past years, various phenomenological waveform models have been developed 
for nonspinning (PhenonA [Uj), aligned-spinning (PhenomB [H], PhenomC [20]) and 
processing (PhenomP [22]) BBHs (these are implemented in the LSC Algorithm Library 
(LAL) [HO])- For the nonspinning binary system we choose to use PhenomA, although 
that is not the latest phenomenological model. Since PhenomA was designed to produce 
only the nonspinning waveforms, the speed of generating nonspinning waveforms using 
PhenomA is faster than the speed of generating those aveforms using other models. 
In addition, PhenomA is expressed by a piecewise function, the inspiral part of the 
IMR waveform can be explicitly defined. Using more recent models might give more 
reliable results, but we expect that the dominant effect on our results is due to the 
difference between the IMR signals and the inspiral-only templates for a given waveform 
model rather than its own reliability. Therefore, the choice in waveform model will not 
significantly affect the results. 

Making use of PN-NR hybrid waveforms, Ajith et al [2] proposed a 
phenomenological waveform model (PhenomA) for nonspinning BBHs defined in the 
Fourier domain by, 


Venom(/) = 


( 1 ) 


The effective amplitude is expressed as 


Aeff = 
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where A is the wave amplitude factor whose value depends on the binary masses and 
five extrinsic parameters determined by the sky location and the binary orientation. 
The effective phase is expressed as 


^'efr(/) = 27r/U + (pc + - + ykr) + ^k) , (3) 

^ k=0 

where tc and fc are the coalescence time and the coalescence phase, rj = mim2/M^ is 
the symmetric mass ratio. In equation ([^, 

C( f, (/ _ + <,2/4 ('I) 

is a Lorentzian function that has a width a, and that is centered around the frequency 
/ring. The normalization constant, w = ^ ^ chosen so as to make Aeff(/) 

continuous across the “transition” frequency /ring- The parameter /merg is the frequency 
at which the power-law changes from /“^/® to /“^/^. The phenomenological parameters 
/merg, /ring, d and /cut are given in terms of M and r\ as 


TrM/merg = Oq 7^^ + 6o + Cq , 
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Figure 1. Left: Fourier domain amplitude of a (normalized) phenomenological 
(PhenomA) waveform starting from 100 Hz for a binary with masses 
Large dots indicate /isco and the phenomenological frequency parameters. Right: 
The detector noise amplitude spectrum y/Sn{f) for Advanced LIGO [31] and the 
characteristic frequencies, /isco (black), /merg (red), /ring (blue) and /cut (green), for 
M = 30, 25, 20,15 and 10 Mq with mi/m 2 = 3 from left. 


TTM/ring =airi‘^ + biri + ci, 

irMa = a2rf -\-b 2 r] + C 2 , 

= as + ^3 ^ + C 3 . (5) 

The coefficients aj,bj,Cj, j = 0...3 and Xk,yk,Zk, k = 0,2,3,4, 6,7 are tabulated in 
table 1 of ng. The left panel of Figure illustrates the Fourier domain amplitude 
spectrum and the characteristic frequencies for a binary with masses (10,1O)M0. In 
the right panel of Figure we show the Advanced LIGO noise curve [31] which we will 
use in this work and compare that with the characteristic frequencies for various binary 
masses. 

In the past studies, the frequency cutoff of inspiral waveforms has been generally 
chosen to be the frequency at the innermost-stable-circular-orbit (isco) of the test mass 
orbiting a Schwarzschild black hole: 

TT Mfisco = 6“^/^. (6) 

As in figure[^ this frequency (/isco) is much smaller than the phenomenological frequency 
/merg- In fhc raugc of M < 30 Mq, the ratios /isco//merg are about 0.54 — 0.64. Thus, 
choosing higher cutoff frequencies than /isco are more efficient to get better overlaps with 
IMR signals. Pan et al m proposed the effective ringdown frequency (/erd = LOTFring, 
where Fring is the ringdown frequency of the effective-one-body waveform model) as the 
frequency cutoff. Taking the detector noise spectrum into account, Boyle et al [32| 
also suggested setting the frequency cutoff to a SNR-weighted average of /isco and /erd • 
They found that such frequency cutoffs are more appropriate than /erd, especially in 
low mass region. In this work, we also explore the effect of the frequency cutoff by 
considering two inspiral template models. We define the IMR waveform model XAilZ, 
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the inspiral waveform models X^erg and Xjgco, respectively as 


XMTl = hphenom(/) 
Xmerg = hp^enom ( /) 
Xjsco = hphenom ( /) 


for / G [/low,/cut], 
for / e [/low , /merg], 
for / G [/low, /isco], 


(7) 


where /low is the low frequency cutoff of the waveforms, that depends on the detector 
sensitivity. 


2.2. Parameter estimation: overlap and confidence interval 


In gravitational wave data analysis, a match between a signal (/*) and a template (ht) 
is expressed by a standard inner product weighted by a detector noise power spectrum 
(5'n) as [HHl IHi] 


{hs\ht) = 4Re r 

J fp 


( 8 ) 


//low Sn{f) 

We adopt the analytic fit to the Advanced LIGO noise curve derived by Ref. [3T] as a 
form, 


Sn{f) = 10-^® 


-4.14 

X 


5x ^ 


+ lll( 


1 — + x ^/2 

1 + xV2 



(9) 


where x — ///o with /o = 215 Hz. This noise curve is described in figure and we 
choose /low = 10 Hz. 

In this work, we describe a single detector analysis and use the normalized 
waveforms, /(/) = h{f)/(h\hy^‘^. Then, since the phase rather than the amplitude is 
the main determining factor in our overlap calculations, we do not take into account the 
five extrinsic parameters (i.e., the luminosity distance of the binary, two angles defining 
the sky position of the binary with respect to the detector, the orbital inclination and the 
wave polarization) in the amplitude factor A. In addition, the inverse Fourier transform 
will compute the overlap for all possible coalescence times at once and by taking the 
absolute value of the complex number we can maximize the overlap over all possible 
coalescence phases (see [35] for more details). In this maximization procedure, we apply 
a nearly continuous time shift by choosing a sufficiently small step size as in [HE]- The 
remaining physical parameters in the wave phase are two mass parameters M and rj 
(the phenomenological parameters are also defined by the mass parameters), and in 
this work we use the chirp mass Me = Mrf^^ instead of M to have better performance 
for our approach. Finally, making use of the normalized signal /^(Ao), where Aq is the 
true value of the signal, and the normalized templates hfiX), where Aj = {Mc,ri}, we 
calculate two-dimensional overlap surface as 


P(A) = max 

tci4^c 


{hfi\o)\hfi\)) 

^{hs{\o)\ks{Xo)){UX)\ht{X)) 


( 10 ) 
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Figure 2. Schematic view showing how to calculate the (1-sigma) confidence 
interval (a). The one-dimensional overlap distribution is calculated by marginalizing 
the original two-dimensional overlap surface. The dotted lines indicate P = 0.99876, 
which corresponds to SNR=20 (see section 2.2 for more details). The one-dimensional 
overlap distribution for r] is given in the inset. We use a binary with masses (15, 3)Mq, 
the exact values of ctMc for this binary are presented in table 1. 


Basically, the above overlap formalism is applied to the context of Bayesian 
parameter estimation. In the high SNR limit, the likelihood (L) can be approximated 
by the overlap surface [37] a^ 

lnL(A) = V(l-i"(A)), (11) 

where p is the SNR calculated by = {hs\hs). From this relation, one might expect 
that the conhdence region of the posterior probability density function is associated 
with a certain region in the overlap surface. 

Following the Fisher matrix (FM) formalism, Baird et al. [38] proposed a method to 
calculate the confidence region from the overlap based on the iso-match contours (IM). 
The connection between the confidence region and the overlap surface can be given by 

P > 1 _ Afc(l -P) 

2p2 

where x|(l — p) is the chi-square value for which there is 1 — p probability of obtaining 
that value or larger and the k denotes the degree of freedom, given by the number of 
parameters. In order to calculate the confidence interval for each parameter, we will 
consider one-dimensional overlap distributions (i.e.. A: = 1), those can be obtained by 
marginalizing the two-dimensional overlap surfaces (this marginalization can be done by 
projecting the two-dimensional overlap surfaces onto each parameters axis). Since the 



:|: Since the likelihood is expressed by a Gaussian distribution for high SNRs, this equation implies that 
the high region of P can be expressed by a quadratic function. 




















IM method follows the FM formalism in determining confidence regions assuming flat 
priors, the validity of this method relies on the Gaussianity of the likelihood (SH]- As 
in equation 0. a Gaussian likelihood corresponds to a quadratic overlap at the region 
given by the SNR, and we found that all overlap surfaces were sufficiently quadratic 
at the region of SNR=20. Therefore, we expect that the IM method is reliable in our 
approach. 

For a SNR of 20, the (1-sigma) confidence interval (cr) of the parameter A is 
determined by the distance between the signal Aq and the template Aj, 


(7 = |hA| = |Ai - Ao 


(13) 


when the parameter value of Aj satisfies P{Xt) — 0.99876. In figure we describe how 
to calculate a from the overlap T’(A). As seen in this figure, the overlap distribution is 
almost exactly quadratic (thus symmetric) in this region. Thus, we only consider the 
one-sided overlaps, where Xt < Aq, in our calculations of a. This choice is also because 
the overlap cannot be calculated in the region beyond the physical boundary rj = 0.25. 
When riQ is very close to 0.25, the overlap surface can be obtained only in the region 
Xt < Aq. 


When applying the IM method to the overlap surface, one should note that a 
one-dimensional confidence interval (a) and a two-dimensional confidence region are 
obtained from different overlap regions. For a given overlap surface with a SNR of 
20, cr calculated by equation (13) correspond to a one-sided width of the original two- 
dimensional overlap contour P = 0.99876 {k = 1) because the marginalised overlap 
distribution can be determined by projecting the two-dimensional overlap surface, while 
the two-dimensional confidence region directly corresponds to the contour P = 0.99714 
[k = 2). Therefore, a is not the same as the one-sided width of the two-dimensional 
confidence region, and generally cr is smaller. 

On the other hand, cr for each parameter can be obtained directly from the original 
two-dimensional overlap surface without marginalizations by means of the Effective 
Fisher matrix method proposed by [361137] . In this approach, the FM can be obtained 
by using an analytic function, that is calculated by fitting directly to the original high 
dimensional overlap surface. O’Shaughnessy et al jlQlllI] showed that the statistical 
uncertainties for mass parameters obtained by this method are in good agreement with 
the results of Bayesian Monte Garlo simulations for a nonspinning and a aligned-spin 
binaries Ha and a precessing binary HU. 


3. Result: statistical uncertainty, fitting factor and systematic bias 

In this section, we assume XAilZ as a complete signal model, and consider XM.1Z 
template model to obtain true statistical uncertainties, and X^erg and Xjsco template 
models to calcultate fitting factors and systematic biases. We examine the efficiencies of 
the inspiral template models for both the detection and the parameter estimation. We 
use Advanced LIGO detector noise power spectrum defined in equation (§. As described 
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Figure 3. Statistical uncertainties in percentage (100 x (J\/X) for nonspinning 
BBHs with a SNR of 20. Red dotted lines indicate the constant chirp masses, 
Me = {4, 5,12}Mq from bottom left. 


above, we choose a high SNR of 20 to ensure that the Gaussian approximation holds for 
IM and FM approaches in calculating statistical uncertainties. However, we note that 
the fitting factor and the statistical bias are independent of SNR. 


3.1. IA4TZ templates: statistical uncertainties 


In order to obtain the statistical uncertainty in parameter estimation, we first consider 
lAiTZ templates. We calculate overlap surfaces for nonspinning BBHs with masses 
mp 2 ^ 3Mq and M < 30Mq, from which the confidence intervals are obtained by 
using the IM method described in figure]^ The percentage uncertainties (100 x cf\/\) 
are summarised in figure In the left panel, the uncertainty contours are overall 
aligned with the constant chirp mass curves (red dotted lines) and slightly misaligned 
in highly asymmetric mass region. The percentage uncertainties for Me range broadly 
from ~ 0.0067% to ~ 0.22% depending on the chip mass of the signal. In the right 
panel, the pattern of contours is overall similar to the case for Me, but the contours are 
more misaligned to the chirp curves in highly asymmetric mass region. The percentage 
uncertainties for rj seem to increase linearly depending on the chirp mass from ~ 0.43% 
to ~ 2.0%. Overall, the accuracy of parameter estimation for gmc/Mc is better than 
that for an/p roughly by 1 — 2 orders of magnitude (c.f. figure 6 in Ref. [12] )• 

In the high SNR limit, on the other hand, the Fisher matrix (FM) can be used to 
estimate the statistical uncertainty (refer to Ref. [Hj and references therein for more 
details). The FM is defined by 


F = 

j- ij 


dh 


dh 


dXj 


A=Ao 


(14) 


where Aq is the true value of the signal. Since Fourier-domain waveform models can 
be expressed by analytic functions of the parameters, the derivatives in this equation 
can be obtained analytically, and only one overlap computation has to be performed 
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Table 1. Statistical uncertainties computed by using the analytic FM and 
the numerical IM method for nonspinning BBHs with a SNR of 20 using the 
phenomenological IMR waveform model. 


mi, m2 

SMq , SMq 

ISMq, 

3Mo 

27Mq 

,3Mo 

15Mo, 

15M© 

Method 

IM 

FM 

IM 

FM 

IM 

FM 

IM 

FM 

(TMc ^ 10 ^ [Mq] 

1.75 

1.72 

23.3 

23.0 

52.3 

51.9 

289 

287 

X 10^ 

10.7 

10.7 

10.2 

10.2 

6.91 

6.98 

49.1 

50.1 


numerically to obtain the FM. Thus, the FM method is computationally much cheaper 
than the IM method. In addition, since the FM formalism is also constructed under the 
assumption of the Gaussian likelihood, we anticipate that that gives the results similar to 
our statistical uncertainties calculated by the IM method. Taking the phenomenological 
IMR waveform model into account, we examined the accuracy of the analytic FM by 
comparing the FM estimations to our results, and found that both results were in very 
good agreement within ~ 2% differences for all binaries considered in this work. In 
table we present the comparison results for several binaries. 


3.2. Xmerg templates: fitting factors and systematic biases 


Next, we consider X^erg as a template model to investigate fitting factors and systematic 
biases. When measuring the match between two different waveform models, the fitting 
factor (FF) or equivalently the mismatch (1 — FF), is widely used [5l [TOl [26l 01]: 


FF = max 


(h,(Ao)|h,(A)) 

V(h,(Ao)|h.(Ao))(ht(A)|hi(A))' 


(15) 


FF is the normalized overlap between a signal waveform hs(Ao) and a set of template 
waveforms ht{tc, fc, \) maximized over t^fc and other parameters Aj. Thus, in this 
work, FF corresponds to the maximum value of the overlaps. 


FF = maxP(Mc, 77 ). (16) 

Men 

In gravitational wave data analysis, FF is used to evaluate the detection efficiency. 
The gravitational wave searches use a bank of template waveforms constructed for 
the corresponding mass range of the systems m 03 0SI 07|. Typically, a template 
bank requires that the mismatch between the templates and the signal does not exceed 
3% [2303 including the effect of the discreteness of the template spacing. In this work, 
we use sufficiently dense spacings in the {Mc-rj) plane in order to avoid this discreteness 
effecf[§} For the two different waveform models, the SNR can be defined by 

p={hs\hty^" = {hs\hsy^"FF. (17) 


§ To obtain FF for one signal, for example, we repeat a grid search near the signal varying the search 
area and the template spacings until we can roughly estimate the size of the overlap contour P = 0.995, 
where P is an overlap weighted by the maximum overlap value in that contour, and finally we find FF 
by performing a 51 x 51 grid search in the region P > 0.995 [49]. 
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Figure 4. Fitting factors for Imerg templates for nonspinning BBHs. The red dotted 
line denotes 77crit (see section 3.2 for more details). 


The detection rate is proportional to the cube of the SNR (thus to the cube of FF). 
Therefore, a FF = 0.97 corresponds to a loss of detection rates of 10%. Inhgure]^ 
we show FFs for nonspinning BBHs for I^erg templates. We hnd that the condition 
FF > 0.97 holds for the binaries in the region M < 24Mq, thus we have Merit ^ 2AMq 
for the detection efficiency. On the other hand, we see that FF depends on the total 
mass overall. This is because the contribution level of the MR phases to the complete 
IMR waveform tends to increase as the mass increases. 

Once FF is calculated by the overlap surface as FF = P(M^®‘^, we define the 
bias of the parameter A by the distance from the true value (Aq) to the recovered value 
(A---), 

hx = A--- - Ao. (18) 

The percentage biases (100 x h\/X) are summarised in figure We find that the bias 
also depends on the total mass overall, those can increase up to ~ 0.3% and ~ 4% for 
Me and 77 , respectively. On the other hand, the recovered values (A^'®'^) are positive for 
all true values (Aq), i.e., b\ > 0 . So, for a given chirp mass, increases with increasing 
770 , and when 77 ^'®'^ is equal to the physical boundary 0.25 ,770 can be equal to some critical 
value (77crit), however 77 ’’®^ cannot exceed 0.25 although 770 increases over rjeut- Thus, in 
the range of 77crit < % < 0.25, we always have 77 ’’®^ = 0.25. In addition, since 77 ’’®*^ is 
fixed at 0.25 in this range, as 770 approaches 0.25, the bias ( 6 ^ = 77 ^'®'^ — Vo = 0.25 — Vo) 
approaches 0 (where also approaches 0). These configurations are well described 
by the contours in the red shaded regions in this figure. Although a value beyond the 
physical boundary implies complex-valued masses, the PN waveforms are well behaved 
for 0 < 77 < 1.0 [32] • For detection purposes, Boyle et al [32] showed that allowing 
such unphysical values, FFs can be significantly improved for the binaries above 30 Mq. 
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Figure 5. Systematic biases in percentage (100 x b\/\) for Imerg templates for 
nonspinning BBHs. The red dotted line denotes r]„it 



Figure 6. Fractional biases {b/a) for Tmerg templates for nonspinning BBHs, where 
a is the true statistical uncertainty given in figure The red dotted line denotes r/crit • 


However, since unphysical values are not allowed in parameter estimation, we only take 
into account physical values for the parameter r|||| 

We have shown that both the statistical uncertainty and the systematic bias increase 
with increasing M. However, as predicted by a simple analytic approach in [SI], the 
bias increases more rapidly than the statistical uncertainty. In parameter estimation, a 
more appropriate quantity can be the ratio of the systematic bias (6) and the statistical 
uncertainty (cr). In figure]^ we present the fractional biases (b/a), where a is the true 
statistical uncertainty obtained by using lAiTZ templates as in hgure We find that 
the fractional bias tends to exceed unity if the total mass is larger than ~ 26Mq. This 
means, the systematic bias (produced by a simplification of the template waveforms by 
taking only the inspiral phase into account from the complete IMR phases) becomes 
larger than the statistical uncertainty (calculated by using the complete IMR template 
waveforms for a SNR of 20). Thus, the critical mass for X^erg templates for the parameter 

II In real Monte Carlo simulations, the range of 7] can be incorporated as a prior, and the posterior 
distribution can be affected by the prior. Thus, when the true value of 7] is close to 0.25, the statistical 
uncertainties tend to be reduced compared to the FM estimations (e.g. see figure 2 of Ref. EOj). 
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Figure 7. Fitting factors for Ijsco templates for nonspinning BBHs. The blue dotted 
line denotes 7?crit- For comparison, we include 77crit (red dotted line) and a FF contour 
(red line) for X^erg templates taken from Figure]^ 


estimation efficiency is obtained as Merit ~ 26Mq. 

3.3. Xisco template: fitting factors and systematic biases 

We also take into account Xisco as a template model, and calculate FFs and biases. In 
this model, we only consider the binaries in the range of M < 24M0. FFs for this 
template model are given in figure We find that the dependence of FF on the total 
mass becomes stronger than the case for X^erg, where the contours are almost exactly 
parallel with the fine of constant M. This is because while /merg is a function of both M 
and rj as in equation ([^ , fisco depends only on M as in equation (§■ Since Xisco waveform 
ends the inspiral phase much more quickly than Xmerg waveform, we also find that FFs 
are significantly reduced compared to those for Xmerg (red line). The valid criterion of 
the template bank to satisfy FF > 0.97 is obtained in the range of M < 15Mq. For Xisco 
templates, therefore, we have Merit 15Mq for the detection efficiency. Ajith [10] also 
found Merit ~ 15Mq for 3.5PN TaylorTl inspiral templates with lAiTZ signals. This 
consistency can be explained by that the inspiral part of the phenomenological model 
has been modeled after the 3.5PN TaylorTl approximant. 

The percentage biases (100 x b\/X) for Xisco templates are summarised in figure 
We find that the biases are several times larger than those for Xmerg templates, thus 
the curve of r^crit (blue dotted line) is placed below that for Xmerg templates (red dotted 
line). A similar work was carried out by Bose et al [2H] using XM.1Z signals and 
3.5PN TaylorTl templates, employing the same numerical formalism in calculating the 
biases as in this work. They took into account M and 77 as the mass parameters, 
and found an increasingly negative bias in M. This is explained by the fact that the 
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Figure 8. Systematic biases in percentage (100 x bx/X) for Xisco templates for 
nonspinning BBHs. The blue dotted line denotes r^crit- For comparison, we include 
T^crit (red dotted line) and a bias contour (red line) for Imerg templates taken from 
Figure 1^ 






Figure 9. Fractional biases {b/a) for Xisco templates for nonspinning BBHs, where a 
is the true statistical uncertainty given in figure The blue dotted line denotes r^crit • 
For comparison, we include r^crit (red dotted line) and a fractional bias contour (red 
line) for Tmerg templates taken from Figure]^ 


templates that give the best fit (i.e., FF) tend to have a smaller M, which tends to 
increase a template’s duration, thereby compensating somewhat its lack of the MR 
phases [28]. While, we use Me and rj parameters, and find that the biases for both 
parameters are increasingly positive. The reason is the same, the binary mass depends 
on both Me and r] (M = and we found that a pair of the recovered values 

(M^^^ > Mco, > r]o) always gives a smaller M than the total mass given by the true 

values (Mco^rjo). For example, a binary with masses (10,5)Mq gives {Mo,Mco,/7o} = 
{15Mo,6.0836Mo,0.2222} and = {14.9685Mo, 6.0845Mo, 0.2268}. 

For a binary with (15,8)Mq, we have {Mq^ Mco^Vo} = {23Mq, 9. 4442Mq, 0.2268} and 
= {22.6119Mo,9.4620Mo,0.2341}. 

On the other hand, two results for br^/r] were given in figure 5 of Ref. [28] in low 
mass region, where one was obtained numerically from the overlap surfaces, and the 
other was calculated by the analytic approximation described in [29]. We find that 
































15 



Figure 10. Comparison between overlaps for three template models. All contours 
correspond to P = 0.99714 (SNR=20), where P is an overlap weighted by FF. 


the two results are quite different, the authors did not give sufficient explanations on 
this apparent disagreement. However, our result for brj/r] can successfully explain the 
features in both results. The contours below r^cut show a pattern similar to that in 
the analytical resuh[^ Above ricrit, the biases are quite small near rjo = 0.25 and do 
not appear to change even as M is increased, this trend is consistent with that in the 
numerical result. 

The fractional biases are summarised in figure and these are also several times 
larger than those for /merg templates. The valid criterion of the template bank to satisfy 
6/cr < 1 (where SNR=20) is obtained in the range of M < 17Mq. Therefore, the 
critical mass for lisco templates for the parameter estimation efficiency is determined by 
Merit ~ 17M0. 


3.4- Statistical uncertainties for inspiral templates and the SNR dependence 


When using inaccurate template waveforms, not only the bias is produced but also the 
statistical uncertainty possibly differs from the true statistical uncertainty. To see this, 
we choose three binaries with masses (6, 3), (10, 5) and (14, 6) Mq, and depict their 
overlap contours for XfATZ (black), Imerg (red) and Xisco (blue) templates together in 
figure [Tot where we define the axes by 5X — b\ so that the three overlaps arrange at (0, 0). 
The contours correspond to P = 0.99714, and following equation (12) this corresponds 
to the two-dimensional confidence region with SNR=20, so they indicate the same level 
of confidence region. Here, we define P as an overlap weighted by FF, 


FF’ 


(19) 


SO that the biased overlap surfaces for X^erg and Xisco templates also have a maximum 
value of 1. We see that the red contours are slightly larger than the black contours, and 
this difference tends to increase as the binary mass increases. However, the red contours 
are nearly comparable in size with the black contours for all models. This implies 


% If we allow to range over the unphysical values, the contours will be smoothly extended to 
the blue shaded region, giving a pattern similar to that in the analytical result in the entire low mass 
region. 
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Figure 11. One-dimensional overlap distributions for the three template models. 
Dotted lines indicate P = 0.99945,0.99876,0.99780 and 0.99505, respectively. These 
lines correspond to the lower boundaries of the overlap curves for calculation of the 
confidence intervals with the SNRs of 30, 20, 15 and 10, respectively (see figure 
for more details). We use a binary with masses (14, 6 )M 0 . Note that for lism 
templates, the blue curves are much wider than the black curves, and the quadraticity 
is significantly broken at SNR = 10. 


that the biased statistical uncertainties are comparable to the true statistical 

uncertainties (a) in our low mass region, but this agreement is not guaranteed at very 
high mass region. On the other hand, the blue and black contours are comparable in 
size for a low mass binary, but the difference between the two rapidly increases with 
increasing M. Therefore, the biased statistical uncertainties for Xjsco can be reliable only 
if the binary masses are sufficiently low. 

So far, we have assumed a fixed SNR of 20. In order to investigate the dependence 
of our results on the SNR, we choose one source binary with masses (14, 6 )M 0 . In 
figure we illustrate one-dimensional overlap distributions calculated by using XAilZ, 
2merg and Xjsco templates, where we also define the x-axis by SX — b\ so that the three 
overlaps arrange at x = 0. The true uncertainty is determined by the black curve, 
and the biased uncertainties are determined by the red and blue curves for X^erg and 
Xjsco templates, respectively. For a given SNR, the corresponding lower boundary of 
the overlap curves is denoted by the dotted line (see figure for more details). In 
the above, we showed that for a given SNR, the biased uncertainties obtained 

by the inspiral templates can be considerably larger than the true uncertainties [a] if 
the binary masses are sufficiently high, and this behavior is much more significant for 
Xjsco templates. Similarly, figure 11 shows that if the SNR is too low, can be 

considerably larger than a even in low mass region, especially for Xjsco templates. One 
can see that the difference between the black and red curves is overall small for both Me 
and rj. However, the blue curves are much wider than the black ones at SNR = 10, and 
the quadraticity is also broken. We found that this discrepancy can be more significant 
for more massive binaries. We therefore conclude that for Xmerg templates are 

overall acceptable in our mass region for any SNRs above ~ 10. For Xjsco templates. 
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SNR SNR 

Figure 12. Comparison between parameter estimation uncertainties (cr) computed 
by the analytic FM and the IM methods with various SNRs (left) and the corresponding 
fractional biases {b/a) for a bias calculated by an IMR signal and inspiral templates 
(right). We assume a binary with masses (16,8)Mq for the signal. 


however, can be acceptable only if the binary mass is sufficiently low and the 

SNR is sufficiently high. 

The FM formalism implies that the true statistical uncertainty (cr) is inversely 
proportional to the SNR. To see this from our overlaps, we calculate the confidence 
intervals using the IM method varying the SNR for a binary with masses (16, 8)Mq. In 
the left panel of figure [T^ we compare the result with the FM result, and find a very 
good agreement between the two methods. On the other hand, FF and systematic bias 
(bx) are independent of the SNR but depend only on the template model (ht). Thus, 
for a bias calculated by an IMR signal and inspiral templates, the fractional bias {b/a) 
increases with increasing SNR, and those are illustrated in the right panel of figure 
for Xisco and Imeig templates. This result indicates that X^erg template model becomes 
much more efficient than Xisco for parameter estimation as the SNR increases. For Xisco 
templates, if the SNR is lower than ~ 12, buc can be smaller than ctmc but the statistical 
uncertainty will not be acceptable because the quadraticity of the biased overlap can be 
significantly broken at such a low SNR as seen in figure 

3.5. Aligned-spinning case 

In this subsection, we show some results for aligned-spinning BBHs. We choose to use 
PhenomC [20] because that is the most recent model for this system. In appendix, we 
briefly describe this model. In order to include the spin effect, PhenomC has the effective 
spin parameter % = (1 + (5)xi/2 + (1 — h)x 2/2 where 5 = (mi — m 2 )/M and Xi = Si/mj, 
Si being the spin angular momentum of the fth BH. In this work, we consider only one 
value {x = 0.5) for the spin parameter, and vary total masses from 12Mq to 24Mq for 
a fixed mass ratio of mi/m 2 = 3. First, we calculate statistical uncertainties for mass 
and spin parameters. Since we have seen that FM is sufficiently accurate in estimating 
the statistical uncertainties for PhenomA, we also adopt FM method for PhenomC. In 
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Figure 13. Statistical uncertainties for nonspinning (Nspin) and aligned-spin (Aspin) 
binaries calculated by using the FM method. We use PhenomA for non spinning and 
PhenomC for aligned-spin models assuming fixed mass ratio {milm 2 = 3) and spin 
(X = 0.5). 




Figure 14. Comparison of statistical uncertainties between IM and FM methods 
for a aligned-spin binary with (15,5) Mq and y = 0.5 using PhenomC model. The 
contours indicate P = 0.99876, and the dotted lines indicate the uncertainties obtained 
from the FM. 


figurewe compare the uncertainties (a) for aligned-spinning binaries with those given 
in figure for nonspinning binaries. We find that (JMc for both systems rapidly increase 
with increasing M, but aMc for aligned-spinning BBHs are a bit larger than those for 
nonspinning BBHs within a factor of 2 . Uncertainties in 77 also slowly increase for both 
systems, but the difference in between the two systems is much more significant 
compared to the case for aMc- This is due to the degeneracy between the symmetric 
mass ratio and the effective spin [38], which is already well known in PN theory [5^153] . 
Thus, if we project the three-dimensional confidence region onto the (77 — x) plane, we 
can see a long ellipse, that has a strong correlation between the two parameters (see, 
figure [^. For reference, we also present where we find that the uncertainties almost 
linearly increase. 

For a sanity check between the IM and the FM methods for the aligned-spinning 
system, we compute the three-dimensional overlap ellipsoid for the binary of (15, 5) Mq 
with X = 0.5, and calculate the statistical uncertainties for the mass and spin 
parameters by using the IM method. In figure 14, we illustrate the overlap contours 
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Figure 15. Confidence regions at SNR=15 and 20 calculated by using inspiral- 
PhenomC templates terminated at /merg (left), and their projections onto Me (middle) 
and rj (right) axes. Total masses are 18, 19, 20 and 2 IM 0 from top to down. Mass ratio 
and spin are fixed to be m.l/m.2 = 3 and x = 0.5. Large dots indicate the recovered 
parameters (M™'^, 


of P = 0.99876 in the (Me — rj) and (% — 17 ) planes, respectively. The doted lines 
indicate the statistical uncertainties obtained from the FM method, i.e., {ctmci ^r], — 

{0.0073,0.0078,0.0213}. We find that the uncertainties for both methods are in good 
agreement for all parameters if we consider the one-sided overlap. However, while the 
two-dimensional overlaps calculated by using PhenomA model are nearly quadratic for 
all SNRs above 10 as described in figure [T^ the three-dimensional overlap are less 
quadratic even at SNR=20. In addition, this behavior can be more pronounced as the 
SNR decrease (e.g., see figure 8 of Ref. [SI!)- Therefore, for the aligned-spinning system, 
one should be careful about the choice of SNR in order to apply the IM and FM methods 
to the parameter estimation. This will be studied in detail in a future work. 

Next, we calculate three-dimensional overlap distributions with IMR PhenomC 
signals and inspiral PhenomC templates. In PhenomC, since the amplitude of inspiral- 
merger phase is defined by one smooth function without a transition frequency, we 
terminate the PhenomC waveforms at /merg defined in equation (j^ to obtain the 
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Figure 16. Fitting factors and biases for nonspinning (Nspin) and aligned-spin 
(Aspin) binaries calculated by varying total masses with fixed mass ratio (mi /m 2 = 3) 
and spin (y = 0.5). Note a sudden jump up of biases between M = 19Mq and 20Mq 
(see figure [T^ . 


inspiral-PhenomC template waveforms. We found that all of the three-dimensional 
confidence regions had the long-thin-curved banana shapes in the (Mc,ri,x) space. We 
give some examples in figure where we show the two-dimensional confidence regions 
by projecting the original three-dimensional overlaps onto the (Mg, 77 ) plane, and the 
one-dimensional overlap functions by projecting the two-dimensional confidence regions 
onto each axis. From these overlap distributions, we find an interesting behavior in the 
variation of biases. The recovered parameter 77 '^®'^) is located at the bottom right- 

hand side of the contours at 18Mq, and that is slightly moved at 19Mq. However, that 
is suddenly moved to the top left-hand side at 20Mq, and again slightly moved at 21Mq. 
This sudden movement is well described in the one-dimensional overlap functions, where 
tends to move from right to left, and 77 ’'®® moves oppositely. It seems that the 
overlap surface becomes bimodal at masses between 19Mq and 20Mq, and then the 
position of maximum overlap suddenly moves from one peak to the other peak. 

By exploring the complicated three-dimensional overlap spaces, we investigate FFs 
and biases for the inspiral PhenomC template^ In figure 16, we compare the results 


The procedure to find the overlap distribution P > 0.995 is the same as in the case for the nonspinning 
system, but we use 31 x 31 x 31 grid points in the final grid search. 
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with those given in figures and for the nonspinning inspiral templates (Xmerg)- In 
the top left panel, we find that FFs for the aligned-spinning system are better than 
those of the nonspinning system for all masses considered here. This improvement is 
just due to the expansion of the parameter space from two to three dimensions. Since 
the masses are strongly correlated with the spin, the addition of the spin dimension 
to the overlap space can increase FF compared to the two-dimensional searches. The 
other three panels show biases for nonspinning and aligned-spinning systems. We find 
that biases are also much larger for the aligned-spinning system. In particular, we see 
that while the biases increase gradually with increasing M for the nonspinning system, 
those for the aligned-spinnig system suddenly jump up at between 19Mq and 20Mq as 
described above. 

4. Summary and future work 

Making use of phenomenological waveform model (PhenomA), we defined the IMR 
model XM.TZ and the two inspiral models X^erg and Xisco, and assumed XAiTZ as our 
complete signal model. We described how to calculate the statistical uncertainties in 
parameter estimation from the overlap surfaces and calculated the uncertainties using 
XMXZ templates for nonspinning BBHs in the range of M < 30 Mq with Advanced 
LIGO sensitivity. We investigated the validity of the inspiral templates in detection and 
parameter estimation, respectively, and provided various crucial values in detail. The 
results of this work can be summarised as 

• For XAiTZ templates, statistical uncertainties overall depend on the chirp 

mass of the system and weakly depend on the mass ratio in highly asymmetric 
mass region. The percentage uncertainties for and r} are 0.0067% — 0.22% and 
0.43%—2.0%, respectively with a SNR of 20, and these results are in good agreement 
with the FM estimates of uncertainty within ~ 2% differences. 

• For Xmerg templates, the valid criterion of the template bank to satisfy FF > 0.97 is 
obtained in the range of M < 24M0. The bias increases with increasing mass more 
rapidly than the statistical uncertainty, that begins to exceed the uncertainty at 
M ~ 26Mq. Thus, Merit ~ 24M0 and ~ 26Mq for the detection and the parameter 
estimation, respectively. 

• For Xisco templates, Merit ~ 15Mq and ~ I 7 M 0 for the detection and the parameter 
estimation, respectively. 

• In our mass region, the biased statistical uncertainties calculated by the 

Xmerg templates are slightly larger than the true uncertainties (a) but overall 
comparable. However, those for Xjsco are acceptable only if the binary masses are 
sufficiently small. The difference between a and tends to increase with 

increasing mass or decreasing SNR. 

We also showed some results for aligned-spinning binaries with fixed mass ratio 
{mi/m 2 = 3) and spin (y = 0.5). We used the waveform model PhenomC, and 
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calculated statistical uncertainties (a\) using the FM method. Aligned-spinning inspiral 
waveforms were obtained by terminating the PhenomC waveforms at /merg- For these 
binaries we fonnd that 

• For IMR PhenomC templates, ctmc for aligned-spinning binaries are a bit larger 
than those for nonspinning binaries within a factor of 2. However, the differences 
in between the two systems are found to be mnch larger compared to the case 
for aMc due to the mass ratio-spin degeneracy. 

• For inspiral-PhenomC templates, three-dimensional confidence regions have long- 
thin-curved banana shapes in the (Mg, 77 ,%) space. FFs for the aligned-spinning 
system can be better than those for the nonspinning system, but biases are much 
larger. In particular, the confidence regions can have bimodal distributions for the 
binaries with masses between 19Mq and 20Mq. 

In this work, we considered limited binary models for the aligned-spinning system, 
so onr results may not be generalised to those binaries with generic masses and spins. 
We will extend our approach to generic aligned-spinning binaries in a future work. 
We showed that the analytic FM method is reliable for estimating the statistical 
nncertainties for the nonspinning system because the overlap surface is nearly qnadratic 
in the (Me, 77 ) plane, However, the three-dimensional overlap was found to be less 
qnadratic even at a high SNR of 20 for our binary model. We will also investigate the 
validity of the FM method for generic aligned-spinning BBHs in detail. The mass-spin 
degeneracy in the aligned-spinning system generally limits our ability to measure the 
individual component masses [SH], but the degeneracy can be broken in the processing 
binaries [HH]- This can also be studied by comparing the FM results for both binary 
systems. 
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Appendix 

The wave amplitude of PhenomC terminates at /cut = 0.15/M, and that is constrncted 
from two parts as 


Aeff — Apu{f)Wj:^ + ApB(f)w^^, 


(A.l) 


where Apm is the premerger amplitude calculated by a PN inspiral amplitude with the 
addition of a higher order frequency term: 


^Pm(/) — ^PN (/) + 


(A.2) 

(A.3) 
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where Q = ttMJ, and Ard is the ringdown amplitude: 

A-rd = ^iC[f, f-KD{a,M),52Q{a)\a)f (A.4) 


where the Lorentzian function is defined by £'(/, /o, d) = d^/[(/ — /o)^ + d/d], and Q is 
the quality factor which depends on the final BH spin a. The two amplitude parts can 
be combined by tanh-window functions: 



1 ± tanh 




(A.5) 


where d = 0.005. The transition frequency /o is determined by /o = 0.98/rd where /rd 
is a ringdown frequency given in terms of M and a. The effective phase is calculated 
by a complete SPA inspiral phasing '/spa, a premerger phasing ■/rm and a ringdown 
phasing '/rd as 


4'efr(/) = /spaM'/i + + /rdwJ, (A. 6 ) 

with fi = 0.1/rd,/2 = /rd using d = 0.005 in the window functions. The premerger 
and ringdown phasing have the forms 

1 


' 0 PM — “(<^i/ + <^2/ ^ + <^3/ + 0^4 + <^5/^^^ + Q^e/)? 

7 ] 

0RD = 01 + 02/, 


(A.7) 

(A. 8 ) 


where the ak coefficients are inspired by the SPA phase, redefined and phenomenolog¬ 
ically fitted to agree with the PN-NR hybrid waveforms, while 0 i ^2 parameters are not 
fitted but obtained from the premerger ansatz by taking the value and slope of the phase 
at the transition point /rd- The coefficients introduced in this model are expressed in 
terms of 77 and y, and those are given in table 2 of j 20 ]. 
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